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ABSTRACT 

Given  a  function  u{x)  which  is  represented  by  its  cell-averages  in  cells  which  are  formed 
by  some  unstructured  grid,  we  show  how  to  decompose  the  function  into  various  scales  of 
variation.  This  is  done  by  considering  a  set  of  nested  grids  in  which  the  given  grid  is  the 
finest,  and  identifying  in  each  locality  the  coarsest  grid  in  the  set  from  which  u{x)  can  be 
recovered  to  a  prescribed  accuracy. 

We  apply  this  multi-resolution  analysis  to  ENO  schemes  in  order  to  reduce  the  number  of 

numerical  flux  computations  which  is  needed  in  order  to  advance  the  solution  by  one  time- 

step.  This  is  accomplished  by  decomposing  the  numerical  solution  at  the  beginning  of  each 

time-step  into  levels  of  resolution,  and  performing  the  computation  in  each  locality  at  the 

appropriate  coarser  grid.  We  present  an  efficient  algorithm  for  implementing  this  program 

in  the  one-dimensional  case;  this  algorithm  can  be  extended  to  the  multi-dimensional  case 

with  cartesian  grids.  - - - -i 
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1.  Introduction 


In  this  paper  we  consider  the  Initial-Boundary  Value  Problem  (IBVP)  for  hyperbolic 
systems  of  conservation  laws  in  s-space  dimensions: 

-h  div/(u)  =  0,  xtD  G  IR^ ,  t>0  (1.1a) 

a(a:,0)  =  2/0(2:),  xeT)  (1-1^) 

with  given  boundary  condition  on  dV,  the  boundary  of  T>.  We  assume  that  the  problem  is 
well-posed  and  denote  its  evolution  operator  by  E{ty,  note  that  it  includes  the  influence  of 
the  boundary  conditions. 

The  computational  domain  D  is  divided  into  cells  {Cj} 

D  =  C,r\Ck^<l>  for  j^k  (1.2a) 

and  we  assume  that  there  is  a  refinement  parameter  h  such  that  the  largest  sphere  contained 
in  each  of  the  cells  is  of  radius  0(h),  and  that  the  ratio  between  the  largest  cell  to  the 
smallest  one  in  the  computational  domain  remains  bounded  under  refinement. 

Let  \Cj\  denote 

\CA  =  f  dx  (1.26) 

JCj 

and  let  uj  denote  the  cell-average  of  u(a:)  over  Cj 

=  MCj)  •  uix)]  (1.3) 

\Cj\  Jc, 

here  A(Cj)  denotes  the  cell-averaging  operator. 

Given  cell-averages  u  =  {uj}  of  u{x)  in  V,  we  denote  by  R(x\u)  an  r-th  order  piecewise- 
polynomial  reconstruction  of  u  from  u,  i.e., 

R{x]u)  =  R,{x-,u)  for  xeC,  (1.4a) 

where  /?,(a;;u)  is  a  polynomial  of  degree  r  —  1.  Expressing  R,{x;u)  as  a  finite  Taylor  series 
around  the  centroid  c,  =  /1(C’,)  •  x 

/{.(x;u)  =  Do  +  '^^  ^(x-c.)'D,  (1.46) 

*=i  i^|=* 

where 
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'  -  A{C,)E{t)R{--v%  =  AiCj)uoix)  (1.5a) 

where  vj  is  an  r-th  order  approximation  to  the  cell-average  of  the  solution  u  at  time 

u"  «  A(Cj)  •  (1.56) 

Due  to  the  diverge  nee- free  form  of  the  PDE  (1.1),  the  scheme  (1.5a)  takes  the  conservation 
form 

"r'  =  /"  /  fW)  ■  Ri-.v'))  ■  Ndsdt  (l.oc) 

|Cj|  v'O  JdCj 

where  dCj  is  the  boundary  of  the  cell  Cj  and  N  is  its  outward  normal.  We  refer  the  reader 
to  [3]  for  details. 

The  purpose  of  this  paper  is  to  present  some  preliminary  results  regarding  the  application 
of  multi-resolution  analysis  to  the  numerical  solution  of  hyperbolic  systems  of  conservation 
laws.  Typically  these  solutions  contain  discontinuities  (shocks,  sliplines,  material  boundaries, 
combustion  fronts)  which  may  move  around  and  also  some  localized  high-frequency  smooth 
behavior  which  is  associated  with  shedding  of  vortices.  In  such  situations  we  have  to  use  a 
fine  grid  in  order  to  resolve  the  details  of  discontinuities  and  vortices,  while  in  other  parts  of 
the  solution  the  use  of  a  coarser  grid  is  adequate. 

The  traditional  solution  to  this  computational  problem  is  to  use  a  nonuniform  adaptive 
grid  where  at  each  time-step  the  discretization  points  are  redistributed  with  the  goal  of  min¬ 
imizing  the  truncation  error  of  the  scheme.  In  the  context  of  the  schemes  (1.5)  this  means 
that  the  cells  Cj  =  are  redefined  at  each  time  and  that  their  size  is  highly  nonuniform: 
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they  are  very  small  around  discontinuities  and  rapid  smooth  variation  of  the  solution,  and  a 
lot  larger  elsewhere.  This  is  certainly  a  reasonable  solution  to  this  computational  problem 
and  there  are  several  such  computer  codes  that  accomplish  this  goal  successfully.  How¬ 
ever,  the  adaptive  grid  approach  is  complicated  and  the  programming  effort  required  for  its 
implementation  is  formidable. 

In  the  following  we  describe  our  concept  of  multi-resolution  analysis  which  is  borrowed 
from  the  realm  of  image-compression  techniques.  There  the  functions  are  defined  on  a 
uniform  fine  grid  and  are  assumed  to  be  over-resolved  in  some  parts  of  it.  The  purpose 
of  the  multi-resolution  analysis  in  image  compression  is  to  determine  appropriate  levels  of 
resolution  for  the  various  parts  of  the  image  in  order  to  eliminate  superfulous  information. 
In  the  context  of  the  numerical  solution  (1.5)  this  means  that  our  computational  grid  (1.2)  is 
rather  uniform  with  respect  to  the  parameter  h  and  is  assumed  to  be  fine  enough  to  capture 
all  the  details  that  we  are  interested  in,  and  that  the  solution  is  in  fact  over-resolved  on 
this  fine  grid  in  large  parts  of  the  computational  domain.  In  order  to  apply  multi-resolution 
analysis  we  construct  an  hierarchy  of  nested  grids  as  follows: 

Given  cells  {Cj}  of  size  hk  and  cell-averages  (1.3)  {u*}  of  u{x)  on  them,  we  define  a 
coarser  grid  {Cj'’*"'}  of  size  hk+i  >  hk  by  joining  some  of  the  cells  in  the  ^'-th  grid  into  a 
single  larger  cell 

C,‘=+'=U,C‘'  (1.6a) 

and  define 

E  IC‘I4  (1.66) 

The  original  grid  is  thus  the  finest  in  the  hierarchy,  and  we  denote  its  cells  by  {C°}  and  their 
size  by  ho  =  h.  To  each  of  these  grids  we  can  apply  the  reconstruction  (1.4)  which  we  denote 
by  R'‘{x;u)  and  also  apply  the  numerical  scheme  (1.5)  in  order  to  advance  its  solution  in 
time. 

What  we  mean  by  multi-resolution  analysis  is  the  assignment  k{i)  such  that 

|/e^<'>(x;u)  - /^"(jr;IZ)|  <  e  for  xeCf  (1.7) 

where  k{t)  is  the  largest  integer  for  which  (1.7)  holds.  In  Section  2  we  shall  describe  a  method 
to  determine  k{i). 

The  general  idea  is  to  use  a  multi-resolution  reconstruction  R{x',u)  which  is  defined  on 
the  finest  grid  by 

R{x\u)  =  R'^^'^x-u)  for  xeC“  (1.8a) 

and  to  consider  the  numerical  scheme 

=  .l(Cf)/-(T)«(.;u").  {LSb) 
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We  observe  that  the  scheme  (1.8)  can  be  thought  of  as  an  adaptive-grid  calculation  in  which 
the  cells  are  a  combination  of  cells  of  various  scales.  We  would  like  to  avoid  this  interpretation 
because  of  its  inherent  complexity,  and  we  look  for  a  simple  implementation. 

In  Section  3  we  develop  such  an  implementation  by  studying  a  related  problem  and  in 
Section  4  we  describe  an  efiicient  multi-resolution  ENO  scheme  for  the  one-dimensional  case. 
In  Section  5  we  summarize  the  present  state  of  affairs  and  point  out  directions  for  future 
development. 

2.  Determination  of  Resolution  Levels 

In  this  section  we  describe  how  to  decompose  a  given  set  of  cell-averages  {uj}  of  a  function 
u(x)  over  cells  {Cj}  of  size  h  into  a  set  of  scales 

h  =  ho  <  hi  <  ■  •  •  <  hi^  (2.1a) 

given  by  the  hierarchy  (1.6).  This  is  done  by  determining  the  larg..st  k  —  k{i)  so  that 

|/f^(x;  IT)  — /f’,(x;  tl)|  <  £  for  xeC,  (2.16) 

for  a  prescribed  level  of  accuracy  c.  Here  /f,(x;  u)  is  the  polynomial  of  degree  r  —  1  in  (1.4) 
and  /?*'(x;I7)  is  the  reconstruction  which  is  associated  with  the  scale  hk. 

First  let  us  describe  the  way  in  which  the  reconstruction  H{x;u)  is  defined  for  the  given 
grid.  Expressing  H,{x;u)  a.s  a  finite  Taylor  series  around  the  centroid  c,  (1.4b)  we  have  to 
find  coefficients  {!)(}. Q  <  K1  <  r  —  1  which  satisfy  the  requirements  (1.4c)  -  (1.4d).  Let  us 
denote  fjy  d  the  vector  of  unknowns 

0<|f|<r-l  (2.2) 

which  are  ordered  in  group.s  of  equal  |fj  with  |f|  =  0, 1,  •  •  ,  r  —  1,  and  denote  by  k  =  ^(r)  the 
number  of  unknowns  in  (2.2).  In  [3]  we  show  how  to  select  a  stencil  of  k  cells  with  indices 
J{i)  (including  ;)  so  that  the  relations 

.4((;’j)  •  /f.(x;  (7)  =  (7j,  jcJ{i)  (2.3a) 

result  in  an  invertible  s}st('m  of  linear  <'qiialions 

Qd  =  a.  (2.36) 


Here  Q  is  the  matrix  with  entri('s 


.4{C’)-(.r-e.)7/^'" 


(2.4a) 
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and  u  is  an  appropriate  ordering  of  the  RHS  of  (2.3a).  Note  that  due  to  the  scaling  in  (2.2), 
the  entries  of  the  matrix  Q  (2.4a)  are  0(1)  under  refinement  and  consequently 

II  <  const,  as  /i 0.  (2-46) 

In  [3]  we  consider  two  cases:  (1)  a  fixed  choice  of  a  centered  stencil  J{i)  which  results  in 
an  upwind  biased  “linear”  scheme  (1.5);  (2)  an  adaptive  choice  of  stencil  J{i)  = 
where  the  stencil  selected  is  the  one  in  which  u[x)  is  the  smoothest  among  several  candidate 
stencils.  This  adaptive  selection  results  in  an  ENO  scheme  (1.5). 

We  turn  now  to  describe  the  decomposition  into  scales  (2.1).  For  A;  =  /\,  —  1,  •  •  •,  we 

check  whether 

max|yl(Cj)  •  R'‘{x;u)  —  Uj\  <  <5(£);  (2-5) 

k{i)  is  the  largest  k  for  which  (2.5)  holds.  Next  we  show  that  (2.1b)  holds  for  an  appropriate 
choice  of  S{e).  Let  us  rewrite  the  polynomial  of  degree  (r  —  1)  /?*'(x;u)  as  a  finite  Taylor 

series  around  the  centroid  c,  as  in  (1.4b)  and  denote  the  vector  ordering  of  its  coefficients 
(2.2)  by  d.  We  denote  by  the  vector  ordering  of 


i‘  =  A(C,).fi‘(iiu),  ]eJ(i) 

(2.6a) 

as  in  (2.3a).  Clearly 

Qd^u" 

(2.66) 

where  Q  is  the  matrix  in 

(2.3b).  It  follows  therefore  that 

Q{d  -d)  =  u-  ii\ 

(2.6a) 

lld-d||<||g-^ll  ||u-  u*^||<||g-^l|.5(£) 

(2.7a) 

and  consequently 

max  |/?^(x;  u)  —  /?,(x;  u)|  <  C\\Q  *||<5(£:) 

x(.C\ 

(2.76) 

where 

r.  F  -  r 

L  —  max  max  . ... 

l<K|<r-l  xeC,  /iKI 

(2.7c) 

and  11  II  is  the  fi-norm.  We  note  that  the  terms  C  and  ||Q~^||  in  (2.7)  depend  on  the  geometry 
of  the  cell  and  of  the  stencil  J(?),  respectively. 

Up  to  this  point  we  have  dealt  with  the  concept  of  multi-resolution  analysis  in  the  context 
of  approximation  of  functions.  We  turn  now  to  consider  its  application  to  the  numerical 
solution  of  the  IBVP  (1.1)  by  the  scheme  ( 1 .5).  Given  {i’"},  the  cell-averages  of  the  numerical 


solution  at  time  over  the  cells  {Cj},  we  apply  the  multi-resolution  analysis  (2.5)  to  find 
k{i)  and  define  the  multi-resolution  reconstruction  R{x',v'^)  (1.8a)  by 

fi(x;u")  = /?*'<‘)(x;u")  for  xeCi  (2.8a) 

and  consider  the  numerical  scheme 

=  A{Cj)E{T)Rix;v^).  (2.8b) 

We  observe  that  due  to  the  well-posedness  of  the  problem 

ll£(r)R(x;  u")  -  E(t)R(x;  u”)||  <  const. ||i?(a:;  v")  -  R(x;  u")ll  (2.8c) 

and  therefore  the  deviation  of  the  values  computed  by  the  multi-resolution  scheme  (2.8) 
from  those  of  the  original  scheme  (1.5)  is  of  order  e.  Roughly  speaking,  all  the  values 
in  (2.8b)  for  which  k(j)  =  ko  can  be  obtained  from  the  calculation  on  the  coarse  grid  ko  by 

uf  ^  =  /l(C,)£(r)R"o(a:;u"),  k{j)  =  ko-  (2.9) 

What  we  need  now  is  an  efficient  algorithm  that  will  enable  us  to  perform  the  computation 
in  (2.8)  using  the  appropriate  coarse  grid  calculations  (2.9)  at  the  computational  cost  of  a 
corresponding  adaptive  grid  implementation  without  its  inherent  complexity.  To  get  ideas 
for  such  an  algorithm  we  consider  a  related  computational  problem  in  the  next  section. 


3.  A  Related  Problem 


In  this  section  we  consider  the  problem  of  computing  discrete  values  of  a  composite 
function  (l){v{x)),  where  v{x)  is  a  function  of  the  type  that  we  encounter  in  the  solution  of 
hyperbolic  conservation  laws  (1.1),  i.e.,  it  is  piecewise  smooth  with  different  scales  of  smooth 
variation.  The  function  </)(t’)  is  a  model  for  the  numerical  flux  in  (1.5c).  We  assume  that 
(^(u)  is  a  smooth  function  of  v  which  is  expensive  to  compute,  e.g., 


/oo 

g(v,x)dx 

■OO 


with  g{v,x)  that  depends  smoothly  on  v. 

In  this  model  problem  we  assume  that  t;(x)  is  defined  in  [0,1]  and  its  discrete  values  are 
given  on  the  uniform  grid 


Cl  —  0,  —  2  •  j,  Vj  —  v(xj). 


(3.2a) 


We  assume  that  the  grid  is  fine  enough  to  resolve  v{x)  to  our  satisfaction,  and  that  in 
fact  v{x)  is  over-resolved  in  some  parts  of  the  grid.  Our  task  is  to  find  an  efficient  algorithm 
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to  obtain  values  (f)^,  0  <  j  <2*^  which  approximate  <f){v{xj))  within  a  specified  tolerance 

of  error  e,  i.e., 

\4>^j  -  <  £■  (3.26) 

We  consider  the  nested  set  of  grids  {6'*^},  0  <  k  <  K,  which  is  defined  by 

=  J>‘ =  {2'‘i}Co'  (3.3a) 

in  which  is  the  finest  (given)  grid  and  is  the  coarsest;  let  6*  denote  the  spacing  of 
points  in  the  grid  G^ 

hk  =  2’‘-’.  (3.36) 

Clearly 

G‘nG*-‘ =G‘  =  {i„}2«j-.  (3.3c) 

G*  '  —  G*  =  ■  (3.3(/) 

The  restriction  of  the  values  {uj}  to  the  grid  G^  constitutes  a  level  of 

resolution  of  v[x)  which  corresponds  to  the  scale  hk- 

We  begin  the  computation  on  the  coarsest  grid  G^  by  calculating  the  values  {4>^]jtji< 
from  the  given  “expensive”  expression,  say  (3.1) 

j^J^-  (3.4) 

Then  for  k  =  A' , .  •  • ,  1  we  use  the  already  computed  values  <f>'^  on  the  grid  G^  to  obtain 
acceptable  approximate  values  of  <^(y(a:))  on  the  finer  grid  as  follows: 

At  all  points  (3.3c)  which  are  common  to  both  grids,  we  retain  the  already  computed 
values 

^  (3.5a) 

The  points  (3.3d)  are  the  ones  added  to  by  splitting  each  of  its  intervals  into  two.  At 
each  such  middle  point,  we  make  a  decision  whether  to  compute  (p  directly,  i.e., 

-  <}>{v2(-i).  {2i  -  (3.56) 

or  to  approximate  p  there  by  interpolation  from  the  grid  G*,  i.e., 

=  /(^2r-i;/),  (2^-  l)cj''-*.  (3.5c) 

We  elect  to  interpolate  (3.5c)  wherever  it  can  be  ensured  that  the  interpolated  v'alue  of  4>  is 
accurate  to  a  prescribed  tolerance 

\<i>\;\-<p{v2e-^)\<ek-^.  (2£-l)cj''-'.  (3.5d) 
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The  basic  idea  behind  this  algorithm  is  that  is  to  be  computed  directly  until  the 

prescribed  level  of  resolution  implied  by  (3.2b)  is  obtained.  Once  the  required  level  of 
resolution  is  achieved  in  a  certain  locality  we  cease  to  use  the  expensive  direct  computation 
and  the  values  of  <f>  are  transferred  to  finer  grids  by  interpolation  (which  is  assumed  to  be 
considerably  less  expensive).  The  basic  assumption  is  that  (f>{v)  is  smooth  and  does  not 
vary  much  in  [miny(x),maxu(x)].  Thus  the  quality  of  approximating  <i>{v{x))  on  by 

interpolation  from  the  coarser  is  determined  up  to  a  scaling  factor  by  the  error  in  doing 
so  for  v{x).  Since  v(x)  is  given  on  all  grids,  our  success  in  interpolating  from  can  be 
easily  measured.  The  criterion  for  acceptable  interpolation  is  formulated  by 

inaxl/(xju:  v^)  -  <  Sk  =>  l/(xj;  -  (l>{vj)\  <  Sk  (3.6a) 

Here  the  “max”  denotes  checking  the  quality  of  the  interpolation  at  the  point  =  0) 

and  possibly  at  neighboring  points  of  the  finest  grid  {I  —  ±1,...).  8k  is  the  corresponding 
tolerance  for  the  v-iriterpolation,  and  it  is  to  be  determined  by  analytical  methods. 

We  recall  that  the  values  of  (f)^  that  are  actually  used  in  the  interpolation  (3.5c)  are  an 
Zk  approximation  to  the  exact  values  thus 

|/(x;<;!>*=)-/(.r;<;^.(i;''-))|<C,.£,,  (3.66) 


and  therefore  in  (3.5d) 

|/(xy;  <^^)  -  <  |/(x,;  -  I{xy,  <^(e'=))|  +  |/(x,;  4>{v'^))  -  <^(i;,)|  <  (1  +  Ck)  •  Sk.  (3.6c) 

Consequently  the  tolerance  levels  satisfy 

fo  =  t  >  ti  >  •  •  •  >  fK  =  0,  (3.7a) 

(1  +  Ck)  ■  i^k  ^  (3.76) 

We  turn  now  to  consider  the  practical  implementation  of  the  multi-resolution  algorithm 
(3.4)  -  (3.6).  The  computational  effort  of  executing  it  is  composed  of  A^(l  —  2“^"^)  checks, 
M  <  N  direct  calculations  of  (p  and  (iV  —  M)  interpolations  oi  (f>\  N  =  2^ .  The  number 
of  checks  is  fixed  and  is  almost  iV,  so  we  want  to  make  the  test  in  (3.6a)  as  simple  as 
possible.  On  the  other  hand  we  wish  to  reduce  A/,  the  number  of  direct  evaluations 
this  calls  for  an  elaborate  test  so  as  to  not  miss  points  in  which  interpolation  is  acceptable. 
These  considerations  lead  to  a  compromise  in  the  level  of  sophistication  to  be  used  in  the 
interpolation  l{x,v^)  and  I{x\(p^).  The  simplest  and  least  expensive  interpolation  is  based 
on  the  use  of  a  fixed  central  stencil.  However,  v{x)  is  a  discontinuous  function  and  using 
a  fixed  central  stc'iicil  means  that  we  shall  miss  all  the  points  near  the  discontinuity  which 
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otherwise  could  be  well-approximated  by  the  more  expensive  ENO  interpolation  that  uses 
an  adoptive  stencil  [1],  [2].  Investing  even  more  in  an  ENO  interpolation  by  using  subcell 
resolution  [4]  v/ill  enable  us  to  get  a  good  approximation  e''en  in  a  cell  which  contains  a 
discontinuity. 

A  reasonable  strategy  is  to  use  a  hierarchy  of  checks:  first  to  try  the  simplest  interpolation 
with  fixed  central  stencil  which  is  a  good  enough  test  for  most  of  the  domain.  If  the  test  in 
(3.6a)  fails  for  the  simple  interpolation  we  may  try  to  get  an  acceptable  approximation  by  an 
ENO  interpolation.  In  cases  where  the  computation  of  <f>{v)  is  so  expensive  that  it  justifies 
an  additional  investment  in  checks,  we  may  try  to  use  also  a  subcell  resolution  approach. 

4.  Multi-Resolution  for  One-Dimensional  Conservation  Laws 

In  this  section  we  present  a  multi-resolution  algorithm  for  the  solution  of  one-dimensional 
hyperbolic  conservation  laws  which  is  motivated  by  the  algorithm  (3.4)  -  (3.6)  for  the  model 
problem. 

We  consider  the  IBVT  for  the  one-dimensional  conservation  law 

««  + /(w)jr  =  0,  0  <  x  <  1,  t>0  (4.1a) 

u(x,0)  =  uq{x)^  0  <  X  <  1  (‘I-16) 

with  appropriate  boundary  conditions  at  x  =  0,.r  =  1. 

We  discretize  functions  in  [0,1]  by  taking  their  cell-averages  on  the  grid  G'°  (3.2a) 

= /l([xj_i,Xj])  •  u(<),  l<i<2'^. 

Given  cell-averages  IT  =  {uj}jfjo  we  consider  the  reconstruction  R{x]u)  which  is  defined 
in  a  piecewise  manner  in  [0,1]  and  satisfies 

R{x;u)  =  u{x)  -|-  0{h''),  wherever  u{x)  is  smooth,  (4.2a) 

A{[xj-i,Xj])  ■  R{x-,u)  =  Uj,  1  <  y  <  2*^.  (4.26) 

We  use  the  numerical  cheme  (1.5)  which  can  be  written  in  this  case  as 

= /l([Xj_i,Xj])  ■  Elr)  • /?(.T;n”),  n°  =  A([xj_i, Xj])  •  Uo(x),  1  <  j  <  2-^,  (4.3n) 

and  has  the  conservation  form 

=  r;  -  Xoif,  -  /,_,),  Ao  =  r/ho  (4.36) 

with  the  numerical  tiux 

7,  =  -  r nEm-ii\.,)dt.  M.3c) 

^  T  Jo  ^ 
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Note  that  the  computation  in  (-1.3)  is  to  be  performed  with  Aq  which  is  restricted  by  the 
CFL  condition,  and  that  /q  and  /2J  include  the  influence  of  the  corresponding  boundary 
conditions  at  x  =  0  and  a'  =  1. 

In  solving  one-dimensional  problems,  as  well  as  in  multi-dimensional  problems  on  carte¬ 
sian  grids,  it  is  convenient  to  work  w’ith  the  primitive  function 

t/(x,0  =  fuU,t)dC  (4.4a) 

Since 

j[U{xj,tn)  -  U{xj.iJn)]  =  v],  (4.46) 

-[U{Xj,tn  +  t)~  U{Xj,tn)\  =  7j.  (4-4c) 

it  is  convenient  to  define  the  reconstruction  Ii{x,v^)  by 

/f(x;c")  =  ^/(x;f/")  (4.5) 

dx 

where  /(.r;6''‘)  is  a  piecewise  interpolation  of  {U(Xj,tn)}  (see  [2])  and  to  compute  the  nu¬ 
merical  fluxes  from  (4.4c). 

We  tur:i  now  to  describe  a  multi-resolution  algorithm  for  carrying  out  the  computation 
of  the  numerical  solution  (1.3).  We  are  given  c'*  =  {r'" on  the  grid  and  our  task 

is  to  compute  all  the  numerical  fluxes  {/^j  on  this  grid  within  a  prescribed  tolerance  for 
error.  Our  basic  assumption  is  that  tlie  solution  is  weil-resolved  on  the  grid  G°,  and  in  fact  it 
may  be  over-  ;esolved  in  some  parts  of  it.  '1  he  computational  task  is  thus  similar  to  the  one 
entailed  in  the  model  problem,  except  that  here  f  ^  is  really  a  functional  of  u(x,t,i)  and  not  a 
function  like  p{v{x))  in  Section  :l.  W(*  use  the  same  set  of  nested  grids  (3.3)  G*,0  <  k  <  K, 
and  initialize  the  algorithm  by  evaluating  in  (4.3c)  on  the  coarsest  grid  G^ .  These  fluxes, 
which  we  denote  by 

7,  =  7(C'),.  (4.6) 

are  actually  an  r-th  order  approximation  to  (4.3c)  w'hich  uses  the  reconstruction  /?(x;u")  of 
the  finest  grid  G^.  Thus  we  use  the  same  numerical  fluxes  that  we  would  use  in  a  straight 
forward  fine  grid  calculation,  except  that  these  fluxes  are  computed  only  in  the  few  locations 
which  correspond  to  the  coarsest  grid  G^' . 

The  algorithm  for  the  numerical  flux  computation  proceeds  in  analogous  way  to  (3.5). 
For  k  —  A'. . .  . ,  1  w(>  use  tin'  already  computed  values  of  the  numerical  flux  /  on  the  grid 
G^  to  obtain  acceptable  values  '  on  th<'  finer  grid  G*“^  as  follows: 

At  all  points  (3.3c)  which  are  common  to  both  grids  we  retain  the  already  computed 
values 

2GJ*-’.  (4.7a) 


10 


At  the  new  middle  points  (3.3d)  we  make  a  decision  whether  to  compute  the  numerical 
flux  directly 

727-1  =7K)2^-i,  (4.76) 

— k 

or  to  approximate  it  by  interpolation  of  /  ,  the  already  computed  values  of  the  numerical 
flux  on  the  grid 

7n\  =  i2£  -  (4.7c) 

The  decision  whether  to  interpolate  or  compute  at  a  certain  locality  depends  on  whether 
we  can  obtain  the  prescribed  accuracy  for  the  numerical  flux  from  the  level  of  resolution 
which  is  offered  by  the  grid  G^.  As  is  indicated  by  (2.8c),  the  error 

7',  -7i=l  [ifim  ■  -  /(£(()  ■  =  £/  (4.80) 

which  is  committed  by  replacing  the  exact  numerical  flux  (4.3c)  by  fj,  the  numerical  flux 
corresponding  to  the  level  of  resolution  which  is  available  in  can  be  bounded  by  the 
corresponding  error  in  reconstruction 


max  \R'‘{x]v’^)  -  R{x;v'')\  =  ef.  (4.86) 

Xj  —  i  + j 

In  the  constant  coefficient  case  f(u)  =  au,  we  get 

|4|<—  r \R\x,-at-,v^)- R{xi-at-,v^)\dt<\a\£f,  (4.9) 

and  similar  estimates  can  be  obtained  in  the  general  nonlinear  case.  As  is  indicated  by  (2.5) 
and  (2.7),  ef'  (4.8b)  can  be  estimated  from 


£“  =  max  |A([a:,_i,x,])  •  i?^(x;  u") 


(4.10a) 


This  check  is  foolproof  but  expensive.  In  our  numerical  computations  we  have  experimented 
with  a  simplified  version  of  (4.10a) 


=max{|A([xj_i,Xj])-/?''(x;u")-i;'’|,  |A([xj, x^+i])  •  i?*(x;  a")  -  1}  (4.106) 

which  takes  into  account  only  the  two  adjacent  fine  grid  intervals  that  determine  the  numer¬ 
ical  flux  at  Xj.  These  numerical  experiments  indicate  that  the  test  (4.10b)  may  be  adequate 
enough. 

From  our  numerical  experiments  we  have  also  learned  that  it  is  better  not  to  interpolate 
the  numerical  fluxes  directly  (4.7c),  but  rather  use  the  following  indirect  way;  Given  we 
can  compute  the  cell-averages  at  time  on  by 


fi^k  7*  >  j,k  pk  •. 

'5 


(4.11a) 
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From  we  get  the  values  of  the  primitive  function  by  using  the  relation  (4.4b). 

Using  interpolation  /(x;  for  the  primitive  function  at  time  fn+j  and  possibly  a  dif¬ 

ferent  interpolation  /(x;U"’'^')  for  the  primitive  function  at  time  we  can  now  obtain  the 
value  of  the  numerical  flu.x  at  the  middle  points  (3.3d)  from  the  relation  (4.4c).  Thus  we 
replace  (4.7c)  by 

757-1  =  {2i-l)ej’^-\  (4.116) 

r 

We  remark  that  unlike  the  model  problem  of  Section  3,  the  computation  of  the  numerical 
flux  involves  propagation  of  discontinuities.  Consequently  a  stencil  of  cells  which  did  not 
contain  discontinuities  of  u[x,  f„)  may  contain  a  discontinuity  of  u(x,  fn+i)  at  its  two  extreme 
cells.  To  overcome  this  problem  we  use  a  central  stencil  for  in  the  test  (4.10b)  and  for  the 
interpolation  I  for  the  primitive  function  at  time  in  order  to  account  for  the  movement 
of  discontinuities  we  take  1  in  (4.11b)  to  be  an  ENO  interpolation  for  the  primitive  function 
at  time 

5.  Summary 

This  work  on  multi-resolution  analysis  has  been  motivated  by  the  current  interest  in 
wavelets.  The  concept  of  wavelets  has  developed  into  a  beautiful  theory  which  is  rich  in 
structure.  The  endowment  of  wavelets  with  many  properties  is  also  their  main  handicap  when 
it  comes  to  practical  computations:  Indeed  wavelets  are  more  local  than  Fourier  analysis, 
but'’  j  ■  •  ‘  r  al'  '  al.  ^!ir  d  sire  to  have  an  orthonormal  basis  with  wavelet  functions 
of  compact  support  results  in  wavelet  functions  which  are  of  fractal  nature.  The  attempt  to 
make  these  wavelet  functions  smoother  results  in  considerably  enlarging  their  support. 

In  the  present  work  we  have  attempted  to  take  from  wavelet  only  its  most  important 
idea  -  the  multi-resolution  analysis  -  and  to  strip  it  of  any  other  structure.  In  doing  so  we 
have  developed  a  multi-resolution  analysis  which  is  devoid  of  beautiful  theory,  but  we  have 
gained  flexibilitv  and  thv  ability  to  adapt  locally.  Given  a  function  u{x)  which  is  represented 
by  cell-averages  on  an  arbitrary  partition  into  cells,  we  have  shown  how  to  decompose  it  into 
various  scales.  This  general  method  provides  a  very  useful  tool  for  removing  superfulous 
information,  especially  when  the  data  is  only  piecewise  smooth,  e.g.,  application  of  this 
methodology  to  image-compression  has  proven  to  be  very  successful. 

VVe  remark  that  this  multi-resolution  analysis  can  be  easily  extended  to  functions  u{x) 
which  are  represented  by  point- values  on  an  unstructured  grid.  The  reason  for  presenting 
only  the  cell-average  version  is  that  we  wish  to  combine  multi-resolution  analysis  with  high 
order  of  accuracy,  and  this  forces  us  to  choose  tlie  scheme  (1.5)  for  accurate  cell-averages. 
We  cannot  use  point-value  schemes  which  are  both  high-order  accurate  and  in  conservation 
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form.  The  need  to  fictitiously  write  in  (4.1a)  as  a  difference  of  a  numerical  flux  results  in 
a  strong  dependence  of  the  numerical  flux  of  the  point-value  scheme  on  the  spacing  of  the 
grid.  Consecpiently  we  cannot  retain  the  already  computed  numerical  fluxes  on  the  grid 
when  refining  to  6'^'“'  (4.7a)  for  point-value  schemes  with  order  of  accuracy  larger  than  2. 

The  multi-resolution  analysis  presented  in  this  paper  has  points  in  common  with  both 
multi-grid  and  adaptive-grid  methods.  It  resembles  multi-grid  in  its  use  of  nested  grids. 
However  multi-grid  is  primarily  an  iteration  scheme  which  is  designed  to  accelerate  con¬ 
vergence,  while  multi-resolution  is  primarily  a  data-compression  scheme  which  is  designed 
to  avoid  superfulous  computations  by  discarding  insignificant  information.  Its  goal  is  the 
same  as  that  of  adaptive-grid  methods  but  it  differs  in  its  concept:  Adaptive-grid  methods 
tell  you  where  to  refine  the  mesh,  while  multi-resolution  analysis  tells  you  where  not.  The 
difference  in  the  praciical  implenientation  of  the  two  concepts  reminds  one  of  the  difference 
between  shock-fitting  and  shock-capturing.  In  adaptive-grid  methods  you  have  to  worry 
about  missing  spontaneously  generated  features.  On  the  other  hand  you  can  do  an  optimal 
job  in  resolving  features  that  you  know  of,  although  at  the  cost  of  complex  programming. 
In  multi-resolution  methods  you  cannot  miss  new  features  because  conceptually  we  always 
compute  on  a  fine  enough  grid.  The  main  issue  is  efficiency:  You  want  to  make  sure  that 
you  do  not  compute  on  too  fine  a  grid,  when  you  do  not  have  to.  An  additional  advantage  is 
the  automatic  way  in  which  the  data  is  handled  in  multi-resolution  analysis  which  translates 
into  simplicity  of  programming. 

The  one-dirnensional  algorithm  of  Section  4  can  be  extended  in  a  straightforward  way  to 
the  multi-dimensional  case  with  cartesian  grids.  Based  on  our  one-dimensional  results  we 
are  led  to  believe  that  realistic  savings  by  factors  between  10  to  100  are  quite  possible  in  the 
three-dimensional  case. 
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